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We use data from the second science run of the LIGO gravitational-wave detectors to search for the grav- 
itational waves from primordial black hole (PBH) binary coalescence with component masses in the range 
0.2-1.0 Mq. The analysis requires a signal to be found in the data from both LIGO observatories, according to 
a set of coincidence criteria. No inspiral signals were found. Assuming a spherical halo with core radius 5 kpc 
extending to 50 kpc containing non-spinning black holes with masses in the range 0.2-1.0 Mq, we place an 
observational upper limit on the rate of PBH coalescence of 63 per year per Milky Way halo (MWH) with 90% 
confidence. 
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Gravitational waves from binary inspiral are among the 
most promising sources for the first generation of gravitational 
wave interferometers. Data from the first and second LIGO 
science runs has been searched for binary neutron star coa- 
lescence with component masses in the range 1-3 Mq 010], 
and a search for binary black holes with component masses 
> 3-M© is underway 0]. Here we consider binaries with 
component masses in the range 0.2-1 Mq. Such binaries 
must contain a pair of black holes in order to be detectable by 
LIGO. Binaries composed of low mass stellar remnants, such 
as white dwarfs, will coalesce before the gravitational waves 
from inspiral reach a high enough frequency to be detected 
by ground based interferometers! 4]. Black holes with masses 
< 1 Mq are assumed to be primordial black holes (PBHs) 
since there is no known mechanism that can produce sub-solar 
mass black holes as a product of stellar evolution. 

There is evidence from gravitational microlensing surveys 
of the Large Magellanic Cloud (LMC) that ~ 20% of the 
Galactic halo is composed of massive compact halo objects 
(MACHOs) with masses 0.15-0.9 Mq 0]. At present the 
explanation of the observed excess of microlensing events is 
controversial. Self lensing of stars in the LMC cannot ac- 
count for all the observed microlensing events 0] and there 
are a number of potential problems with all the events be- 
ing due to white dwarfs in the halo 0] ■ The nature of the 
majority of observed lenses is unknown U and PBHs with 
masses ~ 0.5 Mq have been proposed as possible MACHO 
candidates GjJ QJ. If the MACHOs are PBHs, it will be 
very difficult to determine this using electromagnetic obser- 
vations Jl2ll . If such PBHs formed in the early universe, then 
it has been suggested that some fraction of the PBHs may exist 
in binaries which are coalescing today lll3lfl4ll . If a significant 
fraction of MACHOs are in the form of PBHs, then estimates 
of the rate of PBH binary coalescence suggest that it may be a 
factor of 100 greater than that of binary neutron stars Ill4l.fl^l . 
If this scenario is correct, the PBH binaries are a promising 
source of gravitational waves and the presence of PBHs in the 
halos of galaxies can be confirmed by the detection of their 
coalescence. 

In this paper we report on a search for PBH binaries in data 
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from the second LIGO science run (S2). The data analysis 
techniques used are identical to those used to search for bi- 
nary neutron stars in the S2 data Ql, the only difference being 
in the choice of the search parameters. No inspiral signals 
were found and so we place an upper limit on the rate of PBH 
binary coalescence in the Galactic halo. We compare this ob- 
served rate to that estimated from microlensing observations 
using the model of PBH binary formation proposed in fl5ll . 
Finally we comment on possible future rate limits as the LIGO 
detectors improve towards their design sensitivity. 

Data for the second science run was taken over 59 days 
from February 14 to April 14, 2003. All three LIGO detec- 
tors at the two observatories were operational: a 4 km and a 
2 km interferometer at the LIGO Hanford Observatory (LHO), 
Washington, and a 4 km interferometer at the LIGO Liv- 
ingston Observatory (LLO), Louisiana. These detectors are 
referred to as HI, H2 and LI respectively. During operation, 
feedback to the mirror positions and to the laser frequency 
keeps the optical cavities near resonance, so that interference 
in the light from the two arms recombining at the beam splitter 
is strongly dependent on the difference between the lengths of 
the two arms. A photodiode at the antisymmetric port of the 
detector senses this light, and a digitized signal is recorded at a 
sampling rate of 16384 Hz. This channel can then be searched 
for a gravitational wave signal. More details on the detectors' 
instrumental configuration and performance can be found in 
fill and ifTvtl . In order to avoid the possibility of correlated 
noise sources between the HI and H2 detectors, we only ana- 
lyze data from times when the LI detector is operational. We 
demand that a candidate event be coincident between the LI 
and one or both of the Hanford detectors to reduce the rate of 
background events due to non-astrophysical sources. 

We refer to yfl for a detailed description of the data analysis 
pipeline. Briefly, we used matched filtering with a bank of fil- 
ters constructed using second order restricted post-Newtonian 
templates El El. The bank is designed so that the loss in 
signal-to-noise ratio (SNR) between a putative signal and its 
nearest template is no more than 5% 12011 . Data was filtered 
in 2048 s chunks and times when the SNR p of a template 
exceeded a threshold p > p* were considered candidate trig- 
gers. Noise transients in the data may yield high values of 
SNR, so a time-frequency \ 2 veto 12111 is used to distinguish 
between such events and inspiral signals. The computational 
resources required to perform the search are proportional to 
the number of templates N, which scales as N ~ m mfn > 
where m m j n is the smallest binary component mass in the 
bank. The available resources limited the template bank to 
binaries with component masses above 0.2 M . The num- 
ber of templates fluctuates over the course of the run due to 
changing detector noise, the average value being 14 178 tem- 
plates in the most sensitive detector (LI). The low frequency 
cut off of the search was 100 Hz due to detector noise at lower 
frequencies; the resulting template durations were between 4 
and 56 s. We can determine the time of an inspiral to within 
1 ms, so to be considered coincident, triggers must be ob- 
served within a time window 6t = 1 1 ms between LHO and 
LLO (the light travel time between the observatories is 10 ms) 
and within St = 1 ms between the detectors at LHO. We use 
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the template bank from the most sensitive detector for all three 
detectors and demand that the mass parameters of coincident 
triggers are identical Sm = 0. We demand that triggers in the 
LHO detectors pass an amplitude consistency test. No ampli- 
tude test is applied to triggers from different obervatories as 
the different alignment of LLO and LHO (due to their differ- 
ent latitudes) can occasionally cause large variations in the de- 
tected signal amplitudes for astrophysical signals. Many tem- 
plates may be trigged nearly simultaneously, forming clusters 
of triggers; the trigger with the largest SNR from each cluster 
is chosen for further study; triggers separated by more than 4 s 
are considered unique. 

The sensitivity of the detectors is measured by determining 
the maximum distance to which the detector is sensitive to the 
inspiral of a pair of O.5M PBHs at SNR of 8; that is, the 
distance at which an optimally oriented binary would produce 
an SNR of 8. This distance is refered to as the range of the 
detector. The detectors were at differing stages of progress 
towards design sensitivity during the S2 run and the sensitiv- 
ity of each detector fluctuated over the course of the run in 
response to different noise sources. The average range of the 
detectors during S2 was 704 kpc for LI, 359 kpc for HI and 
239 kpc for H2. As we demand coincidence with the less sen- 
sitive Hanford detectors, the range of the search is limited to 
the neighborhood of the Milky Way, although there are times 
when LI is sensitive to M31. The PBH binary search uses the 
triggered search pipeline described in 0] which takes advan- 
tage of coincidence and the difference in detector sensitivity 
to reduce computational cost. Data from the less sensitive de- 
tectors (HI and H2) is only filtered if a trigger is observed in 
the most sensitive detector (LI). Since we demand Sm = 0, 
the triggered search is functionally equivalent to filtering all 
three detectors with the same template bank and looking for 
coincidence. 

We algorithmic ally select a subset (approximately 10%) 
of the data to be used as playground for tuning the analysis 
pipeline. The playground samples the entire data set so that it 
is representative of the S2 data and allows us to tune our data 
analysis pipeline without introducing statistical bias into the 
upper limit. The goal of tuning the pipeline is to maximize the 
efficiency of the pipeline to detection of gravitational waves 
from binary inspirals without producing an excessive rate of 
spurious candidate events. The false alarm rate of our search 
was set by the available computational resources. If a trigger 
exceeds the SNR threshold, then a x 2 veto must be performed 
at fifteen times the computational cost of a matched filter. This 
limits the SNR threshold to p* = 7 in all three detectors. We 
tune our detection pipeline by attempting to maximize the de- 
tection efficiency of a population of signals which are added 
to the data and then sought. Since we are interested in PBH 
binaries in the halo of the Galaxy, the population we inject is 
distributed as a standard spherical halo with density distribu- 
tion 



where r is the Galactocentric radius and a = 5 kpc is the 
halo core radius. The halo is truncated at r = 50 kpc. The 
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FIG. 1 : 3270 inspiral signals were injected into the data using a uni- 
form distribution in mi and m,2. Each injection found by the pipeline 
is shown with a point and each injection missed is shown by a circle. 
The lines at mi = m,2 = 0.2 M@ show the edge of the template 
bank. When constructing the upper limit on the rate, we only con- 
sider injections that lie inside the region of parameter space covered 
by the template bank, i.e. mi,m,2 > 0.2 Mq. It can be seen that 
the sensitivity of the search is greatly reduced for binaries with a 
component of mass <~ 0.15. 



component masses of the binaries are uniformly distributed 
between 0.1 and 1.0 Mq. Although the template bank is ter- 
minated at a lower mass of 0.2 Mq, we were able to tune the 
search so that it is possible to detect inspirals with component 
masses down to ~ 0.15 Mq, as shown in FigQ] Detection 
efficiencies for mi, m^ > 0.2 Mq were found to be greater 
than 90%, consistent with that expected from consideration 
of the detector sensitivities. We investigated injected signals 
whose masses were inside the template parameter space but 
which were not recovered and found that the loss was due to 
unfavorable alignment of the binary with the detector antenna 
patterns. However, detection efficiency is uniform across the 
region of mass parameter space covered by the template bank. 
For PBH binaries, the x 2 veto was found to be particularly 
powerful. No triggers were observed in the playground data. 

We estimate the background rate for this search by intro- 
ducing an artificial time offset At > 17 s to the triggers com- 
ing from the Livingston detector relative to the Hanford de- 
tectors. We assume the shifted triggers are uncorrected be- 
tween the observatories, however we do not time-shift the two 
Hanford detectors relative to each other as there may be real 
correlations due to environmental disturbances. The triggers 
which emerge at the end of the pipeline are considered a sin- 
gle trial representative of the output from a search if no sig- 
nals are present in the data. By choosing a time shift much 
greater than the light travel time between the observatories, 
we ensure that a true gravitational wave signal will not be co- 
incident in the time shifted data. If the times of background 



5 



10' 



10' 



A S2 final sample 

* Expect background 



10" 



60 65 



75 80 

P 



B5 



HO 



9!) 



FIG. 2: The mean number of triggers per S2 observation time above 
combined SNR p* . The stars represent the expected background 
based on 60 time shift analyses. The shaded envelope indicates the 
standard deviation in the number of events. The triangles show the 
distribution of events from the final S2 sample. 



triggers are uncorrelated between observatories then the back- 
ground rate is entirely due to accidental coincidences which 
can be estimated using the time-shift analysis. A total of 60 
time shifts were analyzed to estimate the background with 
At = ±17 + 10n s, where n = 0, 1, . . . , 29. 

For a coincident trigger, the SNR observed in LI is denoted 
Pl and the coherent SNR observed in HI and H2 is pjj- The 
distribution of background triggers in the (p-^, pn) plane for 
the PBH binary search showed a similar distribution to that 
of the binary neutron star search 0]; the SNR of background 
triggers in the Hanford detectors was typically larger than that 
in the Livingston detector. In order to combine triggers from 
the two detectors, the SNRs of the triggers were combined as 



Pi 



p!/4 
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with any coincident trig gers in the Hanford detectors com- 
bined coherently I2l 12211 . Fig [2] shows the sample mean and 
standard deviation of the expected number of accidental coin- 
cidence events per S2 observation time with combined SNR 
p 2 > pi computed from the 60 time shifts. This can be com- 
pared with the triggers observed by the search to give a visual 
estimate of the significance of the event candidates. 

The pipeline described above was used to analyze the S2 
data. After applying the data quality cuts and the LI instru- 
mental veto described in |2J, and discarding science segments 
with durations less than 2048 s, a total of 375 hours of data 
was searched for binary coalescence. For the upper limit anal- 
ysis, we only considered the non-playground times amounting 
to 341 hours (The extra 2 hours of data in this analysis com- 
pared with fzt is due to a bug fix applied to the data handling 
routines after completion of the analysis in 0|.) The output of 
the pipeline is a list of candidates which are assigned an SNR 
according to Eq. (0. Only three candidates survive in the fi- 



nal sample. All these triggers lie in non-playground data and 
there are no triple coincident triggers. The two loudest coin- 
cident triggers occur when all three detectors were operating, 
but the SNR in HI was too small to cross the threshold in H2, 
so they were accepted as coincident triggers according to our 
pipeline. The third coincident trigger occurred when only the 
LI and H2 detectors were operating. All three triggers had 
values of combined SNR very close to the threshold value of 
p* 2 = 61.25. 

A trigger is elevated to the status of an event candidate if 
the chance occurrence due to noise is small as measured by 
the background estimation. Event candidates are subject to 
further followup investigations beyond the level of the auto- 
mated pipeline to ensure that it is not due to an instrumental 
of environmental disturbance. Fig |2] shows a cumulative his- 
togram of the final coincident triggers versus p 2 overlayed on 
the expected background due to accidental coincidences. The 
final sample of coincident triggers appear consistent with the 
expected background and so we do not believe that they are 
due to gravitational waves. To verify this, further investiga- 
tion of the three surviving triggers was performed. There is 
evidence of transient excess noise in the detectors at the times 
of all three triggers, although the origin of this noise could not 
be conclusively identified. The presence of transient noise, 
the low SNR of the triggers and their consistency with the ex- 
pected background rate due to noise leads us to believe that no 
gravitational wave signals were detected by the search. 

To determine an upper limit on the event rate we use the 
loudest event statistic 12311 which uses the detection efficiency 
at the signal-to-noise ratio of the loudest trigger surviving the 
pipeline to determine an upper limit on the rate. The rate of 
PBH binary inspirals per Milky Way halo (MWH) is 



90% 



2.303 + In P b 

max I 



yr _1 MWH" 



(3) 



with 90% confidence. T is the observation time of the search, 
Nh is the number of Milky Way halos to which the search 
is sensitive at the SNR threshold p* of the loudest trigger 
p max , and P;, is the probability that all background triggers 
have SNR less than p max . This is a frequentist upper limit 
on the rate. For 1Z > TZ 90 %, there is a probability of 90% or 
greater that at least one event would be observed with SNR 
greater than p*. From the background analysis, we estimate 
that Pi, = 0.3 ± 0.1 (statistical error only); however, for this 
analysis we omit the background term by setting Pj, = 1. This 
yields a conservative estimate of the upper limit on the rate. 

During the T = 341 h = 0.0389 yr of data used in our anal- 
ysis, the largest SNR observed was /o^ lax = 67.4. The number 
of Milky Way halos Nh was computed using a Monte Carlo 
simulation in which the data was re-analyzed with simulated 
inspiral signals drawn from the Milky Way halo population 
described by Eq. Q. Although we have some sensitivity to 
the detection of inspirals with components below this mass, 
we restrict our upper limit to the region covered by the tem- 
plate bank 0.2 < mi, m,2 < 1.0 M Q by discarding all injec- 
tions which have a component mass less than 0.2 M Q when 
computing Nh- Fig. [3] shows the value of Nh as a function 
of the loudest event /0^ nax . At p^ aK = 67.4, we find that 
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FIG. 3: The top panel shows the sensitivity in MW halos Nh of the 
search to the target population as a function of the loudest SNR p max . 
The largest SNR observed in this analysis was p^ax = 67.4 mean- 
ing that the search was sensitive to a fraction Nh = 0.95 MWH of 
the halo. The middle panel shows Nh as a function of total mass 
M = mi + 77i2 of the injected signal. The error bars show the statis- 
tical error due to the finite number of injections in the Monte Carlo 
simulation. The lower panel shows Nh as a function of the symmet- 
ric mass ratio -q — mvmilM ■ We can see that the efficiency is a 
weak function of the total mass, as the amplitude of the inspiral sig- 
nal is a function of the total mass. The efficiency of the search does 
not depend strongly upon 77. 
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FIG. 4: The shaded region shows rates excluded at 90% confidence 
by the observational upper limit on PBH binary coalescence pre- 
sented in this paper as a function of total mass M = mi + 7712 
of the binary. The three points show the rates estimated using Eq. 
for halo models S (M — 1.58 Mq), F (M = 0.44 M ) and B 
(M = 1.84 M Q ) of d. 



ponent masses 0.2-1.0 M Q in the Milky Way halo to be 

K 90% = 63 yr -1 MWH _x . 



(4) 



Nh = 0.95 MWH. The various contributions to the error 
in the measured value of detection efficiency are described in 
detail in 0|. In summary, the systematic errors are due to 
uncertainties in the instrumental response, errors in the wave- 
form due to differences between the true inspiral signal, and 
the finite number of injections in the Monte Carlo simulation. 
In this analysis, we neglect errors due to the spatial distribu- 
tion of the PBH binaries as studies show that the upper limit is 
relatively insensitive to the shape of the Milky Way halo. This 
is because the maximum range of all three detectors is greater 
than 50 kpc for PBH masses > 0.2 Mq. The systematic errors 
will affect the rate through the measured SNR of the loudest 
event. We can see that the efficiency of the search depends 
very weakly on the SNR of the loudest event, again due to the 
range of the search compared to the halo radius. The statistical 
errors in the Monte Carlo analysis dominate the errors in Nh- 
The combined error due to waveform mismatch and the cali- 
bration uncertainty is found to be O(10~ 4 ) MWH. The effects 
of spin were ignored both in the population and in the wave- 
forms used to detect insp iral signals. Estimates based on the 
work of Apostolatos|24] suggest that the mismatch between 
the signal from spinning PBHs and our templates will not sig- 
nificantly affect the upper limit. To be conservative, however, 
we place an upper limit only on non-spinning PBHs; we will 
address this issue quantatively in future analysis. Combining 
the errors in quadrature and assuming the downward excur- 
sion of Nh to be conservative, we obtain an observational 
upper limit on the rate of PBH binary coalescence with com- 



By considering numerical simulations of three body PBH 
interactions in the early universe Ioka et. al. 1H5I1 obtain a 
probability distribution for the formation rate and coalescence 
time of PBH binaries. This depends on the PBH mass m, 
which we assume to be the MACHO mass. From this distri- 
bution, we may obtain an estimate of the rate of PBH coales- 
cence at the present time, given by 

R = ixi °" s (^)(^r yr " MwH " ,5) 

where m is the MACHO mass and Ai is the mass of the 
halo in MACHOs, which is obtained from microlensing ob- 
servations. These measured values depend on the halo model 
used in the analysis of the microlensing results l25ll23l . The 
halo model in Eq. corresponds to model S of the MACHO 
Collaboration 1251 . The microlensing observations and PBH 
formation models assume a (^-function mass distribution, as 
does the rate estimate in Eq. Q. We can see from Fig [3] that 
our detection efficiency is not strongly dependent on the ra- 
tio of the binary masses rj, and so we can marginalize over 
this parameter to obtain the rate as a function of total PBH 
mass M, which can be compared with the predicted rates 
from microlensing for different halo models. The analysis of 
5.7 yrs of photometry of 11.9 million stars in the LMC sug- 
gests a MACHO mass of m = 0.79±° 24 and a halo MACHO 
mass M = lOtj x 10 10 M for halo model S 0. Assum- 
ing all the MACHOs are PBHs, we obtain the rate estimate 
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1Z = 1.2 x 10 _2 yr _1 MWET 1 , which is three orders of mag- 
nitude lower than our measured rate. Fig 0] shows a compar- 
ison of the rates predicted using the results of 0] for a stan- 
dard halo (S), a large halo (B) and small halo (F). Models B 
and F are power-law halos 12711 and are discussed in detail in 
l25l Ecill . We note that our upper limit is not strongly depen- 
dent on the halo model as all three halos terminate before the 
sensitivity of our search is significantly decreased. 

Finally we note that the estimated microlensing rate for the 
standard halo is lower than that predicted in tl5ll . due to the 
tighter constraints placed on the MACHO population by the 
additional observation time of |5|]. At design sensitivity initial 
LIGO will be able to see binaries containing 0.5 M@ PBHs 
to 15 Mpc 13, when averaged over antenna pattern and bi- 
nary orientation, suggesting an optimistic rate of several per 
year. The true rate may be much lower, or zero if no PBH 
binaries exist, however the possibility of detection makes it 
worthwhile to extend the inspiral search used for binary neu- 
tron stars into the MACHO mass range. In the absence of 
detection, with 1 yr of data at design sensitivity we should be 
able to place limits on the rate 1Z ~ 10 _3 yr _1 MWH -1 , as- 
suming other galaxies have a similar MACHO halo content 



to our own, and hence significantly constrain the fraction of 
MACHOs that may be PBHs. 
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